function x = g_s(A, b)
    n = length(b);
    A = [A, b];
    for j = 1:n
        % Normalize
        A(j, j:end) = A(j, j:end) / A(j, j);
        
        % Eliminate
        for i = 1:n
            if i ~= j
                A(i, j:end) = A(i, j:end) - A(i, j) * A(j, j:end);
            end
        end
    end
    x = A(:, end);
end